Superfluid-Insulator transition of two-species bosons with spin-orbit coupling 
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Motivated by recent experiments [Y.J. Lin et al, Nature 471, 83 (2011)], we study Mott phases 
and superfluid-insulator (SI) transitions of two-species ultracold bosonic atoms in a two-dimensional 
square optical lattice with nearest neighbor hopping amplitude t in the presence of a spin-orbit cou- 
pling characterized by a tunable strength 7. Using both strong-coupling expansion and Gutzwiller 
mean-field theory, we chart out the phase diagrams of the bosons in the presence of such spin-orbit 
interaction. We compute the momentum distribution of the bosons in the Mott phase near the SI 
transition point and show that it displays precursor peaks whose position in the Brillouin zone can 
be varied by tuning 7. Our analysis of the critical theory of the transition unravels the presence 
of unconventional quantum critical points at t/7 = which are accompanied by emergence of an 
additional gapless mode in the critical region. We also study the superfluid phases of the bosons 
near the SI transition using a Gutzwiller mean-field theory which reveals the existence of a twisted 
superfluid phase with an anisotropic twist angle which depends on 7. Finally, we compute the col- 
lective modes of the bosons and point out the presence of reentrant SI transitions as a function of 
7 for non-zero t. We propose experiments to test our theory. 

PACS numbers: 03.75.Lm, 05.30.Jp, 05.30.Rt 
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I. INTRODUCTION 



Ultracold bosons in optical lattices provide us with a 
wonderful test bed for studying the physics of strongly 
correlated bosons in Mott insulator (MI) and super- 
fluid (SF) phases near the superfluid-insulator (SI) crit- 
ical points. It is well-known that the low-energy 
properties of such bosons can be described by a Bose- 
Hubbard model which captures the essence of the SI 
transition^—. The analysis of such a Bose-Hubbard 
model has been carried out by several authors in the re- 
cent past by using mean-field theory^, quantum monte 
carlo technique^, projection operator method 9 , and 
strong-coupling expansio n 10 ' 11 . The advantage of the 
last method is that it provides a direct access to boson 
Green function in the strongly coupled regime and hence 
to the momentum distribution of the bosons in the MI 
phase near the quantum critical point. In particular, the 
method predicts the occurrence of a precursor peak in the 
momentum distribution of the bosons in the MI phase 
near the SI transition point which has been experimen- 
tally verified^. More recently, several theoretical^ and 
experimental 14 proposals of generating artificial Abelian 
gauge-fields have been put forth. The strong-coupling 
expansion has been also used to describe the SI transi- 
tion of the bosons in the presence of such fields^; such 
studies has also been extended to the case of non- Abelian 
gauge fields 1 ^. Further, the method has also been used 
to study the properties of the bosons in the presence of 
a modulated lattice and it has been shown that such a 
study can reveal the excitation spectrum of the bosons 
both the MI and SF phases near the SI transition point 1? . 

Spin-orbit coupling plays a key role in shaping the low- 
energy properties of several materials including topolog- 
ical insulators which have been a subject of intense re- 
search in recent tiniest. However, the strength of the 



spin-orbit coupling is an intrinsic property of these ma- 
terials and hence not widely tunable. More recently, 
there has been several theoretical proposals of realiza- 
tion of analogous couplings for neutral bosons in a trap 
which has the advantage of generating a tunable spin- 
orbit coupling—. One such proposal has recently been 
realized experimentally^. In the experiment of Ref. I2H 
two suitably detuned Raman lasers was used to gener- 
ate a momentum and spin-dependent coupling between 
the m,F = and mp = — 1 hyperfme states of F = 1 
Rb atoms. These two states acts as two species of the 
bosons and such a coupling is shown to generate a term 
H so = Eek x a y /ki in the Hamiltonian describing these 
atoms. Here Eg = h 2 k 2 /2m is the natural energy unit 
constructed out of the wavelength of Raman lasers ki 
and the mass m of the bosons, and a denotes Pauli 
matrices in the hyperfine space Qmp = 0,-1)) of the 
bosons. We note that such a term is a linear combination 
of the Rashba iJpt ~ (& x k y — a- y k x ) and the Dresselhaus 
Hd ~ ~cr x ky — (jyk x terms. In addition to the spin-orbit 
term, the Raman lasers which are detuned by an energy 
8 from the Raman transition frequency lead to two addi- 
tional terms in the atom Hamiltonian. The first of these 
is directly proportional to the detuning and is given by 
Hd — Sa y /2 while the second term depends on the cou- 
pling strength f2 of the atoms to the lasers: H c — VL<j z /2. 
Together these terms yield an effective Hamiltonian of 
the atoms given by 



H cS = H 2 k 2 I/2m + H so + H d 



(1) 



where I denotes the identity matrix. We note that the 
outset that although the experiments of Ref. [U generates 
H so which is a linear combination of Rashba and Dressel- 
haus terms, there are several theoretical proposals^ - for 
specific generation of either Rashba or Dresselhaus terms 
using Raman lasers. 
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The possibility of realization of spin-orbit coupling for 
neutral bosons has led to several theoretical work on 
the subject - — . Most of these focus on the weak cou- 
pling regime (where the boson interaction can be treated 
perturbatively) and deal with the nature of the possi- 
ble ground states^ 3 -, spin-Hall effect in the presence of 
a shallow tilted lattice and novel spin excitation o 24 ! 27 ' 28 , 
realization of analog of chiral confinement in one- and 
multi-dimensional condensates^, presence of a spin- 
stripe phase 2 ^, dynamics of bosons in the presence of 
spin-orbit coupling using Gross-Pitaevskii equations, na- 
ture of collective excitations 3 ^, and the presence of half- 
quantum vortex excitation a 22 ' 29 of these bosons in the 
presence of the s pin- orbit term in the SF phase. In 
contrast, Refs. l3lU33l focus on the strong-coupling limit 
and derive possible effective spin Hamiltonian to describe 
these phases for fi = 6 = 0. However, the analysis of 
these papers do not provide access to the bosons Green 
functions and do not take into account the effect of finite 
6 and il. One of the central goals of the present work 
constitute obtaining such a Green function in the pres- 
ence of S and f2 and using it for analyzing the critical 
theory of the SI transitions. 

In this work we consider two-species bosons in the pres- 
ence of a spin-orbit coupling term and in a 2D square op- 
tical lattice. The two species of bosons may be thought 
to correspond to two hyperfine states Rb F = 1 atoms. 
In the absence of the spin-orbit coupling and in the pres- 
ence of the lattice, the Hamiltonian for such a two-species 
systems can be written a a 34 i 35 

H = y~"j-(in ia + Uh ia (h ia - l)/2] + XU ^ ^n^a 



(2) 



where b la denotes the bosons annihilation operator on 



the i site, a = 1,2 is the species index, hi, 



b]h„ 



is the boson number operator, {/(AC/) is the intra-(inter- 
) species interaction strength between the bosons, and t a 
(with t\ = t and £2 = ijt) denotes the nearest neigh- 
bor hopping amplitudes. In the presence of the Raman 
lasers inducing a Rashba spin-orbit coupling, the addi- 
tional terms in the boson Hamiltonian are given, in terms 
of a two component boson field ^ = (bn, bi2) T , by 



(3) 



Here the first term represents the lattice analogue of 
the Rashba spin-orbit coupling generated by the Raman 
lasers 3 ^, dij is unit vector along the x — y plane between 
the neighboring sites i and j, O is the species-dependent 
shift in the chemical potential of the bosons, and i5 de- 
notes the detuning as in Eq. [TJ The phase diagram of 
the Hamiltonian given by Eq. [5] has already been studied 



in detail o 34 ' 35 ; the main purpose of this work is to study 
the additional features of the phase diagram due to the 
presence of the terms in Eq. |3] We note here that for 
77 = 1, and (5 = = 0, Hq + H\ is formally equivalent to 
the Hamiltonian studied in Refs. kHT-1331 



The key results that we obtain from such a study are 
the following. First, we chart out the phase diagram of 
the bosons in the Mott phase in the presence of small 
spin-orbit coupling 7 and hopping amplitudes t a . Using 
a strong coupling theory, we also obtain the Green func- 
tion and hence the momentum distribution of the bosons 
in these Mott phases. We find that the momentum dis- 
tribution of the bosons develops precursor peaks near the 
SI transition and show that the position of these peaks 
in the 2D Brillouin zone can be continuously tuned from 
(k x ,k y ) = (0,0) to (k x ,k y ) = (vr/2,±7r/2) by varying 
the relative strengths of the hopping amplitudes and the 
spin-orbit coupling 7. Second, we analyze the SI transi- 
tion and show that the transition, for t a /7 — 0, provides 
an example unconventional quantum critical point in the 
sense that it has an additional mode which is gapped in 
the superfluid phase but becomes gapless at the critical 
point. We note that the presence of such a critical point 
has been theoretically conjectured for hardcore bosons 
with nearest neighbor interactions 3 - 7 -; however, their pres- 
ence has not been demonstrated so far for boson models 
with finite on-site but no nearest-neighbor interaction. 
Third, we chart out the SI phase boundary and study its 
variation as a function of 7 using a Gutzwiller mean-field 
theory and show that the ground state in the presence of 
a finite 7 is a twisted superfluid phase and that the twist 
angle depends on the ratio j/t^. Finally, we compute 
the collective modes of the bosons and demonstrate that 
system undergoes reentrant SI transition which can be 
accessed by varying 7 at a fixed non-zero t. 



The plan of the rest of the work is as follows. In Sec. 
III1 we chart out the Mott phases of the system and com- 
pute the boson Green function and the momentum dis- 
tribution in these phases. This is followed by Sec. IIII1 
where we construct the effective Landau-Ginzburg (LZ) 
functionals for such the SI transitions, and discuss the 
unconventional nature of the critical point for t a /j ~ 0. 
In Sec. IIV1 we use Gutzwiller mean-field theory to chart 
out the superfluid-insulator phase boundary and show 
that the superfluid ground state is a twisted superfluid. 
This is followed by Sec. [V] where we use the LZ func- 
tionals constructed in Sec. |Hl] to compute the collective 
modes of the bosons in the superfluid phases near the SI 
transition. Finally, we present a discussion of the work 
and conclude in Sec. IVI1 
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II. MOTT PHASES AND THE BOSON 
MOMENTUM DISTRIBUTION 

A. Mott phase in the atomic limit 

In this Section, we shall chart out the Mott phases of 
the system in the so-called Mott or atomic limit where 
7 = t a = 0. The Hamiltonian of the system in this limit 
is given by 

-ffMott = 2J[-[M + QSgn(a)}n ia + Un la {n ia - l)/2] 

ia 

+ XJhnhiz + iS (b\ 2 bn - b^h-A , (4) 

i 

where Sgn(a) = ±1 for a = 1,2. Since all the terms in 
the Hamiltonian are on-site, one can choose a Gutzwiller 
like wavefunction \ip M i) = IL fmrn l n i, "2), where n 1(2 ) 
denotes the occupation of bosons of species 1(2) at a 
lattice site r, and compute the energy of the system 
E [{fnm 2 }} = (ipMi\H mo tt\i>Mi)- Further, since the to- 



tal number of particles per site n = n\ + ri2 commutes 
with i?Mott) the Hamiltonian decomposes into different 
sectors labeled by n. Thus, one can separately compute 
and compare the energy functionals E n = E n [{f nin2 }] 
for each n to find the ground state. For n = 0, Eq = 
while for n = 1 the energy functionals reads 

Ei = -(/x + 0)|/i | 2 + (0- M )|/oi| 2 

+^(/o*i/io - fM- (5) 

A similar expression for E2 and E3 can also be written 
down. For n — 2, we find 

/ -2/i + \U iV26 -iy/26 \ 

E 2 = %{ -~iV25 -2fi -2fl + U * 2 

V iy/25 -2fi + 2H + U J 

(6) 

where vf^ = (fn, /20, fo2) T - Similarly for n = 3, one can 
define * 3 = (/12, /21, /03, /3o) T and obtain 



E. 



( -3/i + 2\U + n + U i2S 

-126 -3/i - fi + 2XU + U 

iy/35 

V -iV36 



-iV35 
iV3S 
-3^ + 3n + 3U 

-3fi-3n + 3U ) 



(7) 



The ground state of the system is then determined by 
minimizing E n for a given set of dimensionless parame- 
ters n/U, A, Q/U, andS/U. 

To chart out the phase diagram, we first consider case 
5 = 0. In this case, all the off-diagonal terms in Eq. [5] 
and [7] vanish and one obtains 

£10 = — M - SI, £01 = — M + £11 = _ 2/i + AC/, 
£20 = -2/i - 20 + U, £02 = -2/i + 2fl + U, 
£30 = -3/i - 3ft + 3U, £ 03 = -3/i + 3ft + 3U, 
£21 = -3/i + 2AC7 -0, + U, 

£12 = -3/i + 2AC7 + fi + U. (8) 

The MI phase diagram for 6 = is shown in Fig. Q] 
We note from Eq.[5]that the boundary between MI phase 
(1, 0) and (1, 1) is determined by £10 = £11 leading to the 
condition /it = f2 + AC/. Similarly, the boundary between 
(1,0) and (2,0) phases is determined by the condition 
/x = — CI + U while that between the (1,1) and (2,0) 
phases is given by XU = —2S1 + U. 

For finite S, the energy of different Mott phases are 
determined by Eqs. [3 [6] and [7] Using these equations, 
we find the ground state numerically as function of fi 
and A for several representative values of f2 and 6 as 
shown in Fig. [5] We note that the main effect of S is 



to smoothen out the phase boundary between the phases 
and to realize a MI ground which a linear superposition of 
states with different n\ and n 2 with a fixed n = m + ri2- 
For example, the ground state with n = 1 in Fig. [2] is a 
linear superposition of the states (1,0) and (0,1). The 
overlap of a state (711,712) with the ground state with 
n = m + 772 depends on the precise values of 5 and fi. 



B. Momentum distribution in the MI phase 

In this section, we shall compute the momentum distri- 
bution of the Green function in the Mott phase for which 
n = 1 . The calculations can be generalized to any n in a 
straightforward manner; however this requires handling 
quite complicated algebra which we refrain from in this 
work. 

First, let us consider the Green function of the bosons 
in the MI phase in the atomic limit. For n = 1, the Green 
function is a 2 x 2 matrix given by 



Go(r,r') 



/ (T6 t 1 (r)6i(r')) {Tb\{r)b 2 {r')) . 
\(Tbl(T)h(r')) (TbUr)b 2 (T')) 1 " ' 



To compute the Green function, we first consider the 
eigenenergies of -HMott- These are obtained by diago- 
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FIG. 1: Schematic phase-diagram showing preferred particle 
distribution in /i-A plane for 5 = 0. The pairs of numbers 
(ri.1,712) denotes the particle numbers of the two species at 
each site. 
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FIG. 2: Schematic phase-diagram of the MI phase showing 
in the plane for finite Q, and S. The numbers 1, 2, and 
3 denotes the total boson number n at a given regime. The 
pair of in the box for each panel denote the values of (fl, 6) 
for which the phase diagram has been drawn. 



nalizing E n for the n particle sector; for computing the 
zero-temperature Green function for the n = 1, sector, 
we shall need the expressions of these energies for n = 1 
and n = 2 sectors. For n— 1 , let us denote these energies 
by £ \ and £\ with £ \ < £ 2 . it can be easily seen from 
Eq. [SI that the corresponding eigenstates \£{) and \£ 2 ) 



are related to the states |1,0) and |0, 1) by 

l^i 1 ) N _ 



Ul Vi 

u 2 v 2 



11,0) 

|o,i) 



(10) 



where v\ — u 2 — S/"D, and u% = 
VFTF) /V where V = [5 2 + {Q 2 + VS 2 



i' 2 



W)} 1 ' 2 . Sir 



ilarly for the n = 2 sectors, we denote the cigcnenergies 
and corresponding eigenfunctions of i?Mott by £\ 2 3 and 
|£i 2,3) respectively. From Eq. El we find that the states 
l^iVs) ar e related to |1, 1), |2,0), and |0,2) by 




x\ yi zi 
x 2 V2 z 2 
23 2/3 z 3 




(11) 



where the expressions of Xi, yi and Zi can be found by 
diagonalizing the energy functional E 2 (Eq. [5]). These 
coefficients are found numerically in the present work for 
finite 5. Here we note that Xi, yi and Zi are imaginary 
for i = 2, 3 and real for i = 1. Using these expressions, 
a straightforward calculation following Ref. [ll| yields the 
atomic limit Green functions as 



G° 12 {iw) = 
G° 21 (iw) = 
G° 22 (iw) = 



E 1 



u x vi 



E 1 



U\Vl 



E 1 



tie 

2 



3 



3 



3 



3 



E 1 - 



lie 





- E 1 - 


iw 




1 12 






- E 1 - 


iw 




1 21 






- F l - 


iw 




T 3 

± 22 






- £ 1 - 


iw 



(12) 



where iuj denotes Matsubara frequency and T d ahl for j 
1, 2, 3 and a, b — 1, 2 are given by 



T 1 



T 1 

± 22 



1 21 



2\ui\ 2 \x 2 \ 2 + \p2u\v\x\x 2 + \J2u\v\X\x\ 

+H 2 N 2 , 

I 1 2 1 1 2 + \Fiu\v\x\x-i + \[2u\v\x\x% l 
+2k| 2 |x 3 | 2 , 

T\ 2 = V2\u\\ 2 x* 2 xi + uiWilxij 2 
+V2\v 1 \ 2 x* 1 x 3 + 2v* 1 uix 2 x 3 , (13) 



and T 2 b and T\ are obtained by replacing all XiS in the 
expression of T^ b by yi and z% respectively. Note that, 
when analytically continued to real frequencies using the 
prescription iui — > u + it, G^ (w) is imaginary for i ^ j 
and real for i = j for e = 0. 

The Green functions obtained in Eq. [12] can be easily 
understood as follows. Each term G® b receives contribu- 
tion from a hole branch which corresponds to removal of 
one particle from the Mott state which cost an energy £\ 
in the atomic limit. The other terms represents contribu- 
tion from the different possible particle branches which 
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corresponds to addition of a particle over the ground 
state with n = 1 and cost energies E^ — £\ for j = 1,2,3. 
The poles of the Green functions occur at these particle 
and hole excitation energies. 

To obtain the Green function for finite nearest- 
neighbor terms t a and 7, we follow the procedure intro- 



duced in Ref. [Til First, we define the bosonic fields as 
ip a i(r) = -0 a (ri,r), where a = 1,2, i denote the site in- 
dex of the optical lattice and r is the imaginary time. In 
terms of these fields, the nearest- neighbor hopping and 
spin-orbit coupling terms given by Eqs. [2] and [3] can be 
written as 



(ij) 



1p2j 



(14) 



where we have omitted the r index of the boson fields 
for clarity, /3 = 1/ksT is the inverse temperature and 
ks is the Boltzman constant which will be subsequently 
set to unity. We then write down the coherent state path 
integral for the bosons and decouple the nearest-neighbor 
hopping and spin-orbit coupling terms by two Hubbard- 
Stratonovitch fields Aj(r) = [Ai i (r)A 2 i(r)] T and so that 
the partition function of the bosons can be written as 
(with h = 1) 



rP 
/o 



quartic part of the resultant action is obtained to be 



rP 

Sl = J dT [X) (V'uAV'fa- (A* Via +b.c)) +H } 



E A * A ~ lA ^ 

(ij) 



Mott 



(15) 



Next, we introduce a second Hubbard-Stratonovitch field 
$i(r) = [$ii(r), < &2i('r)] T to decouple the last term in Si 
(Eq.HH). This leads to 



Z = 
S 2 - 



j o dT [j2 (jPiadr^ia + [A?(*<„ - ^ ia ) + h.< 



-H M ott - *i 

(ij) 



(16) 



We note that the field <&i a (r) have exactly the same cor- 
relators as the original boson fields^. With this obser- 
vation, we integrate out the fields Aj a (r) and ipia{i~) to 
obtain an effective action in terms of the field $i a (r). 
The details of the procedure for doing so is elaborated 
in Ref. |Tl]. After some algebra, the quadratic and the 



uj n :a,b 



+A(k)] a6 $ 6 (fc) (17) 
'5>:(r)$ a (r) 







(18) 



where G° denotes the boson Green functions in the 
atomic limit, k = (w„;k), r = (t;x), and A(k) is given 
by 



A(k) 



t±(cosk x + cosky) 7(1 sin 
I 7(— isinfc^ +sinfc y ) t2(cosk x 



- sin fcj,) 

- cos ky ) 



(19) 



In what follows we shall analyze and S^g within 

mean-field theory to obtain the properties of MI and SF 

phases of the bosons. We shall neglect all high order 

terms in the boson effective action which can be shown to 

be irrelevant in the low-energy, low-momentum limits. 

In the MI phase, (<5> a (/c)) = and the boson action, 

(2) 

within mean- field theory, is given by S^g ■ The momen- 
tum distribution of the bosons in the MI phase at zero 
temperature can then obtained from the boson Green 
function Gefr(fe) = [-G a -\uj) + A(k)] _1 as 



i(k) = [ ^TrG cff (fc) 



(20) 



where Tr denotes matrix trace and we have used uj n — > 
lu + it for analytic continuation to real frequencies. 
To evaluate the integral, we note that the integrand 
TrG e ff(A:) is invariant under an unitary transformation; 
consequently n(k) in Eq. [201 can be written as n(k) = 
/" m duTr G d cS (k)/(2n) where G d s (k) is obtained by diag- 
onalizing G s{k) via an unitary transformation and can 
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FIG. 3: The highest hole and the lowest particle bands of the 
bosons for A = 0.4, fi = 0.0117, S = 0.005(7, t = 0.0217, 7 = 
0.049(7, (j, = 0.175(7, and 77 = 0.5. 



be written as 



G d cs (k) = n 



i=l 



[ W - Ei (k)r 



(21) 



where e$(k) are the band energies which are obtained as 
solution of Det[G7 ff 1 (e; k)] = 0, Ai{w, k) are the residue of 
the Green function at the pole u) — £j(k) which has to be 
determined numerically for finite <5, and riband is the total 
number of such bands. The equation for determining 
these bands can be written using Eq. Q21 HH and [T7] as 

{[Go l (e)]n + 2t 1 (cos(k x ) + cos(k v ))} 
x{[G 1 (e)] 2 2 + 2t 2 (cos(k x ) + cos(k y ))} 

= 4sin 2 (fc. y ) 7 2 + [2 7 sin(fc ;c ) + [G^ 1 (e)] 12 \ 2 (22) 

where we have used the fact that [G ( ^" 1 (e)]i2 = [G^ 1 ^)]^ 
and [G _1 (e)]i2 = — i[G ( ^ 1 (e)]i2. We note that for finite 6, 
[G(7 1 (e)]i2 7^ 0. Consequently, Eq. [22] is invariant under 
k y — > — k y but not under k x — > — k x ; thus the energy 
bands satisfy ei(k x ,k y ) = ei(k x ,—k y ) ^ ei(~k x ,k y ). A 
plot of the highest negative and the lowest positive energy 
bands for representative values of parameters is shown in 
Fig. [3J The plot clearly indicates two minima at k y — 
±ky = ±fco and k x = feo- Also, we find that for the above- 
mentioned parameter values n oan d = 8; there are two 
bands with negative and six bands with positive energies. 

To compute n(k), we note that the contribution to 
n(k) comes from the bands for which e(k) < 0. Labeling 
such energy bands as e~(k) and denoting their number 



FIG. 4: The momentum distribution of the bosons in the 
MI phase showing the precursor peaks moving from center of 
the Brillouin zone to (n/2, ±7r/2) with increasing 7/7. The 
plots correspond to fj, = 0.2(7 and have (from top left to bot- 
tom right) (t, 7 ) to be (0.03,0), (0.03,0.02), (0.025,0.025), 
(0.015,0.04), (0.01,0.045), and (0.0,0.048) in units of U. The 
lighter colors indicate larger values of n(k). All other param- 
eter values are same in Fig [3] 



The plot of n(k) is shown in Fig. [4] As expected, we 
find n(k) develops precursor peaks as one approaches the 
SI transition point by increasing t a and 7. This feature of 
n(k) can be easily understood from Eq. [531 and Fig.[3]by 
noting the following points. First, the energy bands are 
independent of k for t a = 7 = (atomic limit) leading 
to a flat n(k). Second, as we approach the SI transi- 
tion, the gap between the highest band with e;(k) < 
and the lowest band with ej(k) > decreases at special 
points (k x ,ztky) in the Brillouin zone. This results in 
peaks of n(k) at these points as we approach the SI tran- 
sition. These peaks are precursors to the SI transition at 
which the bands touch; the position of these precursor 
peaks depend on the ratio 7/tt (for a fixed rf) and can 
be continually tuned from (tt/2,±it/2) to (0,0) as j/ti 
is decreased. Note that since €i(k x ,ky^) — ti(k X: k 
both (k x , ky) and (k x , —ky) correspond to the peak posi 



v ) I 



by n 



band ' 



one obtains the momentum distribution as 



tion; however, since Ci(k x ,k y ) ^ ej(— k x ,k y ) for finite S, 
n(k) need not (and does not) have a peak at (— k x ,ky) 
unless 5 = 0. Numerically, we find = ky = ko for all 
the parameter range we study. Thus our work demon- 
strates that the key effect of the spin-orbit coupling is to 
shift these precursor peaks from (0, 0) to finite momenta 
(fco,±fco) in the Brillouin zone. In the next section, we 
shall investigate the effect of this shift on the SI transition 
point. 



= x^M^kfff a-^-H-^k) 



4=1 



[e,(k)-e-(k)] 



(23) 



where the sum extends over all bands with < 0. 



III. SUPERFLUID-INSULATOR TRANSITION 

In this section, we shall analyze the SI transition for 
two species bosons with spin-orbit coupling. We use the 
strong coupling Green function developed in Sec. Ill Bl to 
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construct an effective low-energy critical theory for the 
transition. This is followed by the analysis of the critical 
theory in Sec. MI Bl 



A. Critical Theory 



with the elements G^i and G\ 2 given by 



-1 



1 



oj + Eq oj + Eq — U oj — E\ 

E = fi + fl, Ei = (i — fl — XU. (24) 



In this section, we analyze the critical theory of the 
superfluid-insulator transition using and (Eqs. 
im and fT8|) derived in Sec. Ill Bl These terms provide the 
microscopic basis for construction of an effective Landau- 
Ginzburg functional for the MI-SF transition. The ana- 
lytical calculations in this section will be carried out for 
5 = for simplicity; however, we shall provide qualitative 
statements for 5^0 case at the end of this section. 

We consider approaching the critical point from the MI 
side. For 5 = 0, the on-site Green function G° is diagonal 



This allows us to write (G°) 1 as a diagonal matrix 



(G 1 



0\-l 



0) as 





F 2 {oj) 



(25) 



, s 2U(oj + E a ) , 

= oj + E ; = ■ ° , F 2 {iv) = u>-E 1 , 



Using Eq. (35] one can write the effective action S^J(5 



S&\8 = 0) 



o,b=l,2 

Fi(uj) + 2t[cos(k x 
27[— i sin(k x ) H 



i + cos(fc y )] 
sm(ky)} 



27 [i sm(k x ) + sin(fcj,)] 
F 2 {oj) + 2t[cos(k x ) + cos(ky)} 



(26) 



where fc = (w,k). Diagonalizing G (fc), we find the two 
eigenvalues to be 



F+(u) +4t[cos(kx) + cos(fc y )] 
±^F 2 _{oj) + 16 7 2 [sin 2 (fc x ) + sin 2 (^)]] (27) 



where F±(oj) = F\(uj) ± F 2 {oj). Thus the quadratic part 
of the effective action of the bosons can be written as 



S^{5 = 0) = - 



a=± 



d 2 kdu> . , . , , , . , 
^$:(&)A a (fc)$ Q (fc)(28) 



where $ 



„+(-) 



(2tt) 3 



+(_) = oi\ ' '$1 + ct^' ^$2 are linear combina- 
tions of the fields $1 and § 2 and af 2 are the components 
of eigenfunctions of G _1 (fc) corresponding to eigenvalues 
of A± given by 



Fx(u) + 2i[cos(fc : 



- cos 



27 [i sin(fec) + 



sin 



iky)}-\±(k) 

n(ky)} 



At the quantum critical point, for 5 = 0, A_(w = 
0, ±k x = ±k y = k ) touches zero which signifies the 
destabilization of the MI phase. The expression for fco 
and the critical values of t and 7 at which this hap- 
pens can be found from the conditions A_(w = 0,±fc x = 



ifc 1t 



and d\-(u = 0,±k x = ±k y = k Q )/dk 



and yields (with F± = F±{u = 0)) 

F + + 8t cos(fc 
F+t 



sin(fco) cos(fco) + 



4 7 2 + 8t 2 



'Fl + 32 7 2 sin 2 (fc ) 
= (30) 



Eqs. l3Q| provide us the position of the critical point and al- 
lows to find "fc/U (t c /U) and ko for any given t/U (i/U), 
fi/U, A, and fl/U. The basic features of the solution 
to Eq. [3D] is as follows. For 7 = 0, the only solution 

of Eq. [30j is fco = and = (\F-\ - F+)/8. As we 
turn of a finite 7, another possible solution emerges at 

k' = arccos(F+<i 2) (7)/(47 2 + 8[^ 2) (7)] 2 ) where t c 2) is the 
solution of F + + 8i c 2) ( 7 ) cos(fc ) = Fl + 32 7 2 sin 2 (fc ). 
Depending on the chosen fi, SI, and A, there is a crit- 
ical value of 7 = 7 0c , at which t^\j) < t\. At this 
value of 7oc, fco shifts to a non-zero value k' . A similar 
behavior may be inferred by choosing a fixed t and by 
slowing increasing 7 to reach the transition. In particu- 
lar we note that in such cases, for t = 0, fco = 7T /2 and 
7c = \jF 1 F 2 /&. A plot of the phase-diagram based on 
Eq. [30] is shown in Fig. [5] The top left panel of Fig. [5] 
shows the MI-SF phase diagram in the fj, — 7 plane for 
specific t while the top right panel exhibit the phase di- 
agram in the t — fx plane for specific 7. These plots are 
qualitatively similar to their mean-field counterparts in 
Fig. [7] and exhibit reentrant SI transition as a function 
of 7 for any non-zero t. The bottom panels of Fig. [S] 



8 




0.02 0.03 0.04 0.05 0.0 0.01 .,..0.02 0.03 
y/U t/U 



FIG. 5: The phase boundary in the fi — 7 (left panels) and 
jj,—t plane (right panels) for r\ = 1, and Q — 0.01(7 as obtained 
from the strong-coupling analysis. The top panels have 8 = 
while the bottom panels have 8 = 0.005(7. The values of t for 
the left panels corresponding to different lines are t/U = 0.0 
(black solid), t/U = 0.01 (blue dashed), 0.02 (red dotted). 
For the right panels, 'y/U = (black solid line), 0.03 (blue 
dashed line) and 0.04 (red dotted line). 



1.6 n 




t»/u t/u 

FIG. 6: Left panel: Plot of fco against t a (fj.) for 'y/U = 
0.01(pink dash dot dotted), 0.02 (yellow dash dotted ), 0.03 
(blue dotted), 0.0333 (green dashed) and 0.04(black solid) 
line. Right panel: Plot of the line with z = 1 quantum phase 
transition in the 7 — t plane as obtained from solution of 
K Q = O.For all plots 8 = 0, A = 0.4, r) = 1, and fi = 0.01(7. 



quadratic action 



shows the phase diagrams for finite 5 (computed numer- 
ically starting from the expression of Go for finite 8 in 
Eq. and using the method outlined in this section) 
which are seen to be qualitatively similar to their (5 = 
counterparts. The left panel of Fig. [6] shows a plot of fco 
as a function of t c (fi) for several representative values of 
7 with A = 0.4, and Q = 0.01(7. We find that for small 
7, there is a finite range of t c for which the transition 
takes place at fco = 0. The width of this region shrinks 
with increasing 7 and beyond a critical 7 = j' c , the tran- 
sition always takes place with finite fco- For A = 0.4 and 
il = 0.01(7, we find j' c ~ 0.033J7 as can be seen from the 
left panel of Fig. [51 



The critical theory for the MI-SF transition can now be 
constructed in terms of the low-energy excitations around 
lu = and k = fco which can be described by a set of 
bosonic fields ifi(k) around each of these minimum. For 
n such minima at fco 
the field $_ — 3>_(w : 



: fc where i = l..n, one expresses 
0, k = fc ) + (fii(k) and obtain the 



s 2 cr 

Ki 



J SSE^W (Ko^ + K^-v^^k) 



( 27 r)3 



<9A_ 



duj 
Id 2 A 



uj— 0,k— kl 



F-{0)F!_{0) 



F+(0) +8tcos(fc ) 



2 dw 2 



uj— 0,k— kX 



1 r 

4 



^+(0) 



F_ 2 (0) + ^"(0)F_(0) 
F+(0) +8tcos(k l ) 

1 9 2 A- = _ 2 

2 9(fc ) 2 w=o,k=fcj 

8 7 2 cos(2fcj,) 



[Fi(0)F_(0)] 



[F+(0) + 8tcos(fc^)] 3 
t cos(fc ) 
64 7 4 sin 2 (2fc^) 



F+(0) + 8tcos(fc ) [F+(0) + 8tcos(fc?,)] 3 



(31) 



where ' denotes differentiation with respect to to. From 
Eq. [3l] we find that the critical theory has dynamical 
critical exponent z = 2 except for special points at which 
Ko = leading to z = 1. In usual MI-SF transition this 
point appears to be at the tip of the MI lobe. Here we 
find a line of such z = 1 transitions in the t — 7 plane 
as shown in the right panel of Fig. [6] for representative 
values 6 = 0, A = 0.4, 77 = 1, and O = 0.01(7. 

The structure of the quadratic part of the critical ac- 
tion found in Eq. [31] remains qualitatively similar for 
(5^0 except for two differences. The first difference 
in the effective action comes from the fact that the num- 
ber of minima is halved due to the lifting of k x — > —k x 
symmetry as discussed in Sec. Ill Bl while the second dif- 
ference stems from the fact that v x ^ v y for 5^0 leading 
to an anisotropic dispersion of the critical theory. Con- 
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sequently, the critical action S2 now has the form 



s: 



ci;S=£0 



(2.)3 £ 



-v 2 k 2 



v y y ) 



<Pi{k) 



(32) 



The positions of the z = 1 line in the [i — 7 plane also 
changes. However, rest of the features remain the same. 
In Sec. IIIIB1 we shall analyze the critical theory in de- 
tails and show that the MI-SF transition at t/7 = is 
unconventional in the sense that it is accompanied by the 
emergence of an additional gapless mode at criticality. 



B. Analysis of the critical theory 

Having established the analytical form for S" for 5 = 
0, we shall now analyze the effective Landau-Ginzburg 
theory for the transition. The quadratic part of the ef- 
fective action remains the same as in Eq. 1311 Our anal- 
ysis shall hold for 5 ^ as well; in this case Eg . I3T1 shall 
be substituted by Eq. [33] In this section, we shall not 
bother with microscopic calculation; instead we shall an- 
alyze the critical theory from the symmetry perspective 
as done, for example, in Ref. l37l . for small i/7 where the 
minima of G^ 1 {k) occurs at non-zero = (ko,±ko). In 
the presence of two such minima, the bosonic field can 
be written as 



Mr.t) =¥>i(r,t)e i: 



^ 2 (r,t)e l 



(33) 



Substitution of Eq. [33J in Eq. [55] and subsequent expan- 
sion in uj and k (around kj) leads to Eq. 13 II with n = 2. 

To obtain the quartic action, one can in principle sub- 
stitute Eq. [33J in Eq. [T5] average over the fast oscillating 
components involving various powers of cos(k^ • r) and 
sin(kjj : • r) which appears in the expression of S4, and ob- 
tain an effective critical action in terms of </?i,2- However, 
such an averaging proves to be tricky when ir/(k^a) do 
not turn out to be a small integer since one may have 
to sum over an arbitrary large number of lattice sites 
for achieving a proper averaging 1 ^. Also, for irrational 

, such an averaging procedure is ill-defined. For our 
case, since k^ is a continuous function of t/7, we adopt a 
symmetry-based general method for deriving the fourth 
order term in the action. 

The symmetry based derivation of the effective action 
relies on the fact that an effective low-energy Landau- 
Ginzburg action describing a phase transition must be 
invariant under the projective symmetry group (PSG) 
transformation of its underlying lattice 37 . The elements 
of PSG for a square lattice are translation by a lattice 
vector along x and y (T x and T y ), rotation about the z 
axis by 7r/2 (R w / 2 ), and reflection about x and y axes (P x 
and P y ). Following the method derived in Ref. [33 and 
using Eq. 1331 we find that under these transformation the 



bosonic field ip(r, t) transforms as 
T 3 



ipi -> e^Vi, 



V2 -> e- ik ° a ip 2) 



T„ : pi -> e^V, V2 - 
P x : if! -> ip 2 , (p 2 -> ipi 

Py-.lfl^ Ifl, lfi2 -> fl 



e lfeo V 2 , 
-> e 2lfeo > 2 



(34) 

To find the fourth order effective action consistent with 
Eq. [34] we first consider the case t — for which ko = 
7r/2. In this case, the most general form of the quartic 
action is 



d 2 rdt 



(M 



1^2 1 2 ) +770(^1^2 +h.c) 2 
(35) 

where 770 is a constant whose value will be determined 
later. Redefining the fields £x( 2 ) = (jfii + (— )f2)/V^, one 
gets 



ncr 9 

64 " 2 



d 2 rdt 



(l&l a + l6lT + »(lfi| a -|6l 2 ) a 



(36) 

For r]o > 0, the ground state of Sf and S 2 thus cor- 
respond to condensation of both the fields: £1 = 
and £2 = £2 ex P(*Mo)- However, the relative phase Ho 
between these two fields is not fixed by the S" . In- 
deed, if we construct the eighth order term S$ in the 
effective action, it will have a PSG allowed term S$ — 
A' J d 2 rdt{Q&+h.c) 4 which will fix /_t = rmr/2 for A' < 
and fj,o — (m + l/2)n/2 for A' > where m is an inte- 
ger. Thus the effective phase mode characterized by the 
fluctuation of the relative phase /Jo is massive in the SF 
phase but is expected to become gapless when A' — > 
due to irrelevance of S$ at criticality. Consequently, we 
expect all transitions with 770 > to have an additional 
gapless mode in the critical region. To compute the value 
of 770, we note that since n/(koa) = 2, it is possible to 
compute the effective action S" by direct substitution of 
Eq. [33J in Eq. [T51 followed by averaging over fast oscil- 
lating terms as shown in Ref. I la , This procedure yields 
Eq. [35] with 770 = 1. Thus we find that for t/7 = 0, the 
two-species bosons with spin-orbit coupling undergoes an 
unconventional phase transition at 7 = j c which are ac- 
companied by emergence of an additional gapless mode 
at the transition 3 - 7 -. We note that this also implies that 
the vortices corresponding to any one of these fields ip% 
or if 2 will have a fractional vorticity in the sense that a 
boson wavefunction would pick up a phase 7r when moved 
around such a vortex 3 ^. However, generating such vor- 
tices experimentally in present systems may turn out to 
be difficult. 

For t/7 7^ where ko ^ 7r/2, wc find that the only 
form of the effective action which is consistent with the 
PSG transformation has the form 



s 4 cr = 



d 2 rdt 



(M 2 + |H 2 ) 2 +77^i| 2 |H 



(37) 
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The value of jf is difficult to determine for arbitrary fco; 
however, for certain values of fco which satisfies ko/(ira) £ 
Z, one can determine r]' . In all such case we find rj' > 0. 
This indicates that for all t/y, only one of the fields ip\ 
or (f2 condenses. Thus the MI-SF critical points for such 
finite t/y are conventional. 



IV. MEAN-FIELD ANALYSIS 

In this section, we use a Gutzwiller wavefunction to 
obtain the mean-field SI phase boundary for the system. 
The Gutzwiller variational wavefunction which we shall 
use is given by 

= a i |l,0) < + 6 4 |0,l) < + c i |l,l) < + di|2,0) < 
+e i |Q,2) i + / i |0,0) i , 

i* G > - n>*>- ( 3g ) 

i 

Note that for the purpose of charting out the phase dia- 
gram and for describing the SF phase near the SI tran- 
sition point, it is not necessary to incorporate the higher 
number states since we expect these states to have very 
small overlap with the ground state of the system as can 
be checked by explicit numerical calculation. The varia- 
tional energy of the system can be easily computed using 
Eqs. HH and [35] and yields 

E = (* G |(Ho+Hi|* G > =E + E 1 + E 2 
Eo = ^-(M + ^)|a 4 | 2 -(M-^)l^| 2 + (A-2 A1 )|c l | 2 

i 

+ (1 - 2/j - 2fl)\d l \ 2 + (1 - 2(j, + 2fl)\ ei \ 2 
+i5[bid* — dib* + V2ad* — \[2c*di 
+V2c*e l -V2c l e*}, 
Ei = -ti^(A* 1 A jl ,+? 7 A* 2 A j2 )+h.c. 



Eo — — 



m 

7 [£(A* 1 A j - !d2 -A* 2 A^ 1 
+*X)(AS l A i ,a + AJ,A J> i) 



h.c. 



(39) 



where denotes sum over both x and y neighbors of 

site i while Ylla ) denotes sum over x{y) neighboring 
sites of i, and the order parameter A ia — (b ia ) can be 
expressed in terms of the Gutzwiller wavefunction coeffi- 
cients as 



An = f*a t + b* Cl + V2a*di 
A l2 = /*&i + aJci+V2&Jei 



(40) 



The phase diagram obtained by numerical minimiza- 
tion of Eq. [39] is shown in Fig. [Jj We note that for these 
bosons, SI transition can be induced either by tuning y 
or t. We first consider the case of S = 0, 77 = 0.5 and 



0.4-1 



0.3 



30.2 



0.1 



0.0 



• ■ 
• ■ 
• ■ 



0.4 -I 



0.3 



a0.2 



0.1 



0.0 



0.02 
y/U 



-, 0.0 



0.04 0.0 



0.02 
t/U 



0.04 



FIG. 7: Left panel: The MI-SF phase boundary from mean- 
field theory for different representative values of t/U = 0.0 
(black squares), t/U — 0.02 (blue circles) in the fj, — 7 plane. 
Right panel: The phase boundary in the fj,— t plane for y/U = 
0(black squares) and y/U = 0.02 (blue circles). Here we have 
taken 77 = 0.5, 5 = 0.0 and fi = 0.01(7. 



fl = 0.01C/. The MI phase for the parameter values is 
characterized by n\ = 1 and n 2 = 0. The MI-SF phase 
diagram, in the /x — 7 plane, is shown in the left panel of 
Fig. [7] for representative values of t/U — 0, 0.02. Here, 
we find that for all values of /i, the transition always takes 
place into a SF phase with (61), (62) 7^ 0. Following the 
nomenclature of Ref. HH, we term this SF phase as 2-SF. 
We also note that for any finite t, the bosons display reen- 
trant SI transition with variation of strength of 7. The 
MI-SF phase diagram in the fj,— t plane for representative 
values of 7 = 0, 0.02 is shown in the right panel of Fig. [7] 
Here for 7 =^ 0, we find that the transition always takes 
into a 2-SF phase. In contrast, for 7 = 0, a small region in 
the phase diagram near fi = exhibit 1-SF superfluidity 
for which (61) ^ and (b 2 ) — 0. The phase diagram with 
small non-zero 5 turns out to be qualitatively similar. 

The most striking point about the superfluid phase 
into which the transition takes place becomes evident on 
examining the values of A^ a for the ground state configu- 
ration in the SF phase. We find that although the ampli- 
tudes of the superfluid order parameters remain homoge- 
neous, their phases vary with positions for finite y/t; in 
other words, the superfluid ground state realized is an ex- 
ample of a twisted superfluid phase 38 . We also note that 
the relative phases between the x and the y neighboring 
links are different leading to an anisotropic twist. To ob- 
tain an qualitative understanding of the role of spin-orbit 
coupling in the realization of such a twisted superfluid 
phase, we note that these phases contribute to the energy 
of the system through the terms E\ and E 2 in Eq. [SHI 
Taking cue from the numerical result that the magnitude 
of the order parameters remain constant in the ground 
state configuration, we now write Ai a = Ag a exp(i(fii a ). 
In what follows, we choose the phase of the order param- 
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eter on the i and the neighboring sites as 



9n 



0. 



l 3a I 



?Q 2 = Pj a , 



(41) 



where the subscript a takes values x and y. Using this, 
one can write Ei and E 2 in terms of the relative phases 
between the x and y neighbors 



E 2 

LL *01 



-^2{[cos{a jx ) + cos(a jy )] 
+?7K 2 [cos(/3 Jx - O ) + cos(/3j v - O )] j 

(•i.) 

(ijy) 



(42) 



where k = A02/A 



01- 



Next, we define relative phases living on x and y links 
of the 2D square lattice as 



$42 = <f>j x 2-(f>0, $42 = <j>j y 2 ~ 00 

$^3 = 0^1" 00- (43) 

In terms of these phases, Eq. 221 can be recast as 



~3> 



43 



1(2 



-35, 7t _*,7lt _-2L,2L n 57t 
" 4 — 2 7 4, 2 u, 4 



_Jt 7Jt i n n '5j[ 7t 
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FIG. 8: Left panel: Plot of the relative phases <I>£ X H 1,2,3 on 
the links of the square lattice as a function of 7/i for 77 = 0.5, 
Q. = 0.01£/. Right panel: Distribution of the phases 0i a for 
the superfiuid order parameter for 7/i > 1. 



d(Ei + E-2)/d<5>i xW 1(2,3) = which yields 

sin($£ y i) - — cos($ £x2 - $^3 + $41) = 

r)K 2 sin($4 2 ) - y { sin($^ 2 - $£. E 3 + $41) + cos( 

$42 - $£ x 3 + $41) - COs($ £y2 - $£ x 2 + $4s)} = 

r]K 2 sin($^ 2 ) - -j sin($^ 2 - $42 + $43) = 
sin(5v x3 ) + sin($ £x2 - $43 + $41) + cos($42 

-$4,3 + $41) ~ COs($£ y 2 - $£ x 2 + $43) = (45) 

In general, these equations need to be solved numerically 
and we have not found analytic solutions for them for 
arbitrary values of k and 7/i. However, in the special 
case 7/t ^> 1 , we find that these equations admit an easy 
solution 



Ei 



E 2 
tA 2 01 



~ {[ cos (^,i) +cos($^i)] 

4,4 

+77K 2 [cos($^ x2 ) + cos($^ 2 )]} 

[cos($ £x3 ) - cos($^ x2 - $<? x3 + $41) 



1 K 

T 



4,4 



sin($ £x2 - $£ x3 + $41) 
sin($ £y2 - $42 + $43) 



(44) 



From Eq. [44] we clearly see that unless 7/i is small, the 
minimal energy configuration correspond to non-zero but 
uniform values relative phases over x and y links. Note 
that the precise numerical values of these phases depend 
on k and hence requires input from numerical minimiza- 
tion of Eq. [351 However, once we know the value of k, 
we find that the relative phases for the minimum energy 
are the solutions of the coupled transcendental equations 



$41 = -$41 = 7r / 4 > $4,3 = 7T 
$£ x 2 = -$42 = 3tt/4 



(46) 



The corresponding phase distribution of the superfiuid 
order parameter 0i a = (0jij0i2) is shown in the right 
panel of Fig. [5] For all values of 7/t and k, we find 
$41(2) = -$41(2)- Also, for tj < 1 (which implies 
K < 1), we find that $^ x 2 have a discontinuous jump 
to finite value around 7 = 0. The occurrence of this 
can be easily understood a competition between second 
and the third set of terms ( those proportional to rjn 2 and 
"fK/t respectively) in Eq.UU A plot of the relative phases 
on the x and y links is shown in the left panel of Fig. [8] 
as a function 7/i. We find that the relative phases take 
finite value for non-zero 7 and approaches those given 
by Eq. 3^1 with increasing 7/i thus leading to the real- 
ization of a twisted superfiuid ground state. We have 
checked that the value of the relative phases obtained 
from minimization of Eq. 3?] agrees to those computed 
from minimization of Eq. 1391 
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V. COLLECTIVE MODES 

In this section, we use the critical theory developed in 
Sec. IIII Al to obtain the collective modes in the superfluid 
phase near the critical point. We first consider the case 
S = for which one can obtain straightforward analytical 
expressions for these modes. We first consider t = 0. We 
begin with the quadratic and quartic parts of the boson 
action in the SF phase near the critical point which are 
given by 



S 



So 



S2 ~\~ 5*4 

d 2 kdw 



/d kdio v—v „ / 



i=l,2 

-v 2 \kf + \ ro \) 6 



(47) 



||^[(i6i 2 + i6i 2 ) 2 + %(iai 2 -i6i 2 ) 2 ] 



In the SF phase, b oth the fields condense with amplitudes 
1 60 1 = y\ra \/2g for i = 1,2. To obtain the collective 
modes, we therefore expand the fields £j = + 6£i, 
where <5£; represents small amplitudes fluctuating fields 
which describes the collective modes of the condensate. 
Using Eq. 271 we obtain an effective quadratic action for 
5\i- It turns out that for 770 = 1, the quadratic actions 
for <5£i and 5£,2 reduces to block-diagonal form which can 
be written as 



S' = 



A = 



d 2 kduj x - 



A5( W ,k) + M ^ 

Jai D (-w,-k) 



r 



(48) 



where % = (<5£*,<5&) , and D (w,k) = ~{K Y u? + K uj- 
v 2 |k| 2 + |ro|)/2. The collective modes corresponding to 
the field £j can then be obtained from the condition 
DetA = and yields, 



Wl(2)(k) 



'±a k + ^ai + A{Ai-\r^)Kl 



V 2Kf 
= 2K 1 A^ + K 2 A k = « 2 |k 



M (49) 



Each of these two modes are doubly degenerate. It is easy 
to see from Eq. that CJ2 (k) are gapped while oj\ (k) 
is gapless with tj\ ~ |k| 2 at small |k| for Kq ^ and 
lui = \\s\v/*/K~i for Kq — 0. The mass of the gapped 
mode can be read off from Eq. 05] and are given by 



= J{2\r Q \K x + K*)/K* 



(50) 



Note that in this case, there is one gapless and one 
gapped mode and each of these are doubly degenerate. 
This leads to two gapless modes in the SF phase which is 
a consequence of condensation of both £1 and £2 at the 
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FIG. 9: Left upper panel: Plot of the mass of the gapped 
mode for t = and 8 = 0.005J7 as a function of (7/70 — 
1). Upper right panel: Plot of the velocities of the gapless 
modes (v x and v y ) as a function of (7/7= — 1). Lower left 
(right) panel: Mass of the non-condensed modes (Velocity of 
the gapless mode) for t/7 = 0.005, K = and 5 = 0.005(7. 



Next we consider the case t/7 7^ 0. 
begin with the action 



In this case, we 



S 

So 



So ~\~ Sa 



f d 2 rdt Vt (K 0J + K x u? - t> 2 |k| 2 + \r \) ^ 

J i=l,2 
Sl = \J d 2 rdt [(|^|2 + |^|2)2 + ^|^|2 b2 |2]_ 



(51) 



To obtain the collective modes, we note that th e field <p\ 
condenses with an amplitude |y>oi| = VrolTff- We then 
expand the fields <pi = ipoi + Sipi and cf2 = S<fi2 and 
obtain the effective quadratic action for the field 5tpi t 2- 
It turns out that these actions decouple. The effective 
action for Sipi turns out to be analogous to Eq. 05] and 
yields a gapless and a gapped mode with uj — wi(2)(k). 
The effective action for 5tf2 is given by 
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S = 

a" = 



d 2 kduj i // . 

Do(«,k) + (2 + ^)|r |/2 

D (-£J,-k) + (2 + ^)|r |/2 



(52) 



r 



where \£ 2 = [¥>2(w, k), (/j^— u;, — k)] T . The collective 
modes obtained using Eq. [52] are given by 



W 3 (4)W = (53) 



The masses of these modes are given by 



W 3 (4) 



-(+)K + y/K* +4^1^1(1 + ^) 



(54) 



Thus in this case, we have one gapless and three gapped 
mode. We note that since Kg, K\ and v can be computed 
from microscopic parameters of the theory, our analysis 
provides a way of obtaining the velocities and masses 
of the gapped and the gapless collective modes directly 
from the parameters of the microscopic Hamiltonian of 
the bosons. 

The inclusion of finite 5 changes this picture in two 
essential ways. First, it lifts the degeneracy between some 
of the modes. Second, it makes the dispersion anisotropic 
since in the presence of a finite 5, v x and v y are not 
identical. A plot of the masses of the gapped and velocity 
of the gapless modes for a finite but small S = 0.005J7 is 
shown in Fig. [9] In accordance with the expectation, we 
find that the velocities of the gapless modes are different. 



elusions of our work are the following. First we have 
shown, via explicit calculation of the boson momentum 
distribution function, that the SI transition is accompa- 
nied by precursor peaks in the MI phases near the tran- 
sition and that the position of these peaks can be tuned 
by tuning the strength of the spin-orbit coupling. We 
note that this feature of our theory can be directly ver- 
ified experimentally by routine momentum distribution 
measurement s 1 ' 1 2 . Second, we have analyzed the MI-SF 
phase boundary and have shown the existence of reen- 
trant SI transitions at fixed t and 77 with variation of 
7. This feature can also be detected experimentally by 
momentum distribution measurements. Third, we have 
shown that for t/"f = 0, the SI transition is unconven- 
tional in the sense that it is accompanied by emergence 
of a gapless mode in the critical region. Fourth, we have 
computed the collective modes in the SF phase near the 
transition. We have presented analytical formulae for 
the gapless and the gapped mode and have provided ex- 
plicit expression for their masses and velocities in terms 
of microscopic parameters of theory. These predictions 
can be verified by routine spectroscopy measurements on 
these systems^. Finally, our mean-field study has re- 
vealed the presence of a twisted superfluid ground state 
in these systems with an anisotropic twist angle whose 
magnitude depend on 7/t. 



VI. DISCUSSION 



In this work, we have studied the SI transition of two- 
species bosons with spin-orbit coupling. The main con- 
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